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Formulations and algorithms for the multiple depot, fuel-constrained, 

multiple vehicle routing problem 
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Abstract —We consider a multiple depot, multiple vehicle 
routing problem with fuel coustraiuts. We are glveu a set of 
targets, a set of depots aud a set of homogeneous vehicles, oue 
for each depot. The depots are also allowed to act as refueliug 
statious. The vehicles are allowed to refuel at auy depot, aud 
our objective is to determiue a route for each vehicle with a 
minimum total cost such that each target is visited at least 
ouce by some vehicle, and the vehicles never run out fuel 
as it traverses its route. We refer this problem as Multiple 
Depot, Fuel-Coustralued, Multiple Vehicle Routing Problem 
(FCMVRP). This paper presents four new mixed integer linear 
programming formulations to compute an optimal solution for 
the problem. Extensive computational results for a large set of 
iustances are also preseuted. 

Index Terms —fuel coustraiuts; greeu vehicle routiug; electric 
vehicles; mlxed-iuteger lluear programmiug; uumauued vehicle 
routiug 

1. Introduction 

In this paper, we extend the classic multiple depot, mul¬ 
tiple vehicle routing problem (MDMVRP) to include fuel 
constraints for the vehicles. We are given sets of targets, a set 
of depots, and a set of vehicles, with each vehicle initially 
stationed at a distinct depot. The depots also perform the 
role of refueling stations, and it is reasonable to assume 
that whenever a vehicle visits a depot, it refuels to its 
full capacity. Given this, the objective of FCMVRP is to 
determine a route for each vehicle starting and ending at 
its corresponding depot such that (i) each target is visited at 
least once by some vehicle, (ii) no vehicle runs out of fuel as 
it traverses its path, and (iii) the sum total cost of the routes 
for the vehicles is minimum. Some of the applications for the 
FCMVRP are path-planning for Unmanned Aerial Vehicles 
(UAVs) [1], [2], [3], routing for electric vehicles based on 
the locations of recharging stations [4], [5], and routing for 
green vehicles [6]. Some of these application domains are 
elaborated on the following sections. 

A. Path-planning for UAVs 

Small UAVs are being used routinely in military applica¬ 
tions such as border patrol, reconnaissance, and surveillance 
expeditions, and civilian applications like remote sensing, 
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traffic monitoring, and weather and hurricane monitoring [7], 
[8], [9]. Even though there are several advantages due to 
small platforms for UAVs, there are resource constraints due 
to their size and limited payload. It may not be possible 
for a small UAV to complete a surveillance mission before 
refueling at one of the depots due to the fuel constraints. 
For example, consider a typical surveillance mission with 
multiple vehicles each starting at a depot and together are 
required to monitor a set of targets. To complete this mission, 
the vehicles might have to start at their respective depot, 
then visit a subset of targets and reach one of the depots 
for refueling before starting a new route for the rest of the 
targets. This can be modeled as a FCMVRP with the depots 
acting as refueling stations. 

B. Routing problem for green and electric vehicles 

Green vehicle routing problem is a variant of the Vehicle 
Routing Problem (VRP) and was introduced by authors in 
[6] to account for the challenges associated with operating 
a fleet of alternate-fuel vehicles (AFVs). The US transporta¬ 
tion sector accounts for 28% of national greenhouse gas 
emissions [10]. Several efforts over many decades focusing 
towards the introduction of cleaner fuels (e.g. ultra low 
sulfur diesel) and efficient engine technologies have lead to 
reduced emissions and greater mileage per gallon of fuel 
used. Government organizations, municipalities, and private 
companies are converting their fleet of vehicles to AFVs 
either voluntarily to alleviate the environmental impact of 
fossil based fuels or to meet environmental regulations. For 
instance, FedEx, in its overseas operations, employs AEVs 
that run on biodiesel, liquid natural gas, or compressed 
natural gas. A major challenge that hinders the increase 
in usage of AEVs is the number of alternate-fuel stations 
available for refueling. The ECMVRP is a natural problem 
that arises in such a scenario. An algorithm to compute an 
optimal solution to the ECMVRP would generate low cost 
routes for the vehicles, while respecting their fuel constraints. 

Increasing concerns about climate changes and rising 
green house gas emissions drive the research in sustainable 
and energy efficient mobility. One such example is the 
introduction of electrically-powered vehicles. One of the 
main operational challenges for electric vehicles in transport 
applications is their limited range and the availability of 
recharging stations. The number of electric stations in the 
US is a mere 9,571 with a total of 24,631 charging outlets 
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Fig. 1. Electric station locations in Texas, USA til] 


[11]. Fig. 1 shows a map with the locations of the electric 
stations in Texas, USA; observe that the distribution of the 
electric stations is very sparse except in the four major 
cities Dallas, Houston, Austin, and San Antonio. Successful 
adaption of electric vehicles will strongly depend on the 
methods alleviating the range and recharging limitations. 
If we consider the range and the recharging stations for 
the electric vehicles as analogues to the fuel capacity and 
refueling stations of vehicles that run on fossil-based or 
alternate fuels respectively, then the problem of electric 
vehicle routing subject to the range constraints and lim¬ 
ited availability of electric stations can be modeled as an 
FCMVRP. Clearly, any feasible solution to the FCMVRP 
can be used to implement a feasible route for an electric 
vehicle. 

II. Related work 

The FCMVRP is NP-hard because it contains the VRP 
as a special case. The existing literature on the FCMVRP is 
quite scarce. The multiple depot, single vehicle variant of the 
FCMVRP was first introduced by authors in [12]. When the 
travel costs are symmetric and satisfy the triangle inequality, 
authors in [12] provide an approximation algorithm for this 
variant. They assume that the minimum fuel required to 
travel from any target to its nearest depot is at most equal to 
Fal2 units, where a is a constant in the interval [0,1) and 
F is the fuel capacity of the vehicle. This is a reasonable 
assumption as, in any case, one cannot have a feasible tour if 
there is a target that cannot be visited from any of the depots. 
Using these assumptions, Khuller et al. [12] present a (3(1-1- 
q;))/(2(1 — a)) approximation algorithm for the problem. 
Authors in [1] formulate this multiple depot single vehicle 
variant as a mixed-integer linear program and present A:-opt 
based exchange heuristics to obtain feasible solutions within 
7% of the optimal, on an average. Later, Sundar et al. [2] 
extend the approximation algorithm in [12] to the asymmetric 
case and also present heuristics to solve the asymmetric 
version of this variant. Furthermore, variable neighborhood 
search heuristics for FCMVRP with heterogeneous vehicles, 
i.e., vehicles with different fuel capacities, are presented by 
Levy et al. [3]. More recently, an approximation algorithm 
and heuristics are developed for the FCMVRP by the authors 


in [13]. 

Variants of the classic VRP that are closely related to the 
FCMVRP include the distance constrained VRP [14], [15], 
[16], [17], [18], the orienteering problem [19], [20], and the 
capacitated version of the arc routing problem [21], [22]. The 
distance constrained VRP is a special case of the FCMVRP 
with a single vehicle and single depot that can be considered 
as a fuel station. The FCMVRP is also quite different and 
more general compared to orienteering problem where one is 
interested in maximizing the number of targets visited by the 
vehicle subject to its fuel constraints. Lastly, the arc routing 
problem is a single depot VRP given a set of intermediate 
facilities, and the vehicle has to cover a subset of edges along 
which targets are present. The vehicle is required to collect 
goods from the targets as it traverses the given set of edges 
and unloads the goods at the intermediate facilities. The goal 
of this problem is to find a tour of minimum length that starts 
and ends at the depot such that the vehicle visits the given 
subset of edges, and the total amount of goods carried by 
the vehicle does not exceed the capacity of the vehicle along 
the tour. One of the key differences between the arc routing 
problem and the FCMVRP is that there is no requirement 
that any subset of edges must be visited in the FCMVRP. 

The aim of this paper is to introduce and compare four 
different formulations for the FCMVRP and present branch- 
and-cut algorithms for the formulations. The first two for¬ 
mulations are arc-based, and the rest are node-based formu¬ 
lations that use the Miller-Tucker-Zemlin (MTZ) constraints 
[23]. The major contributions of this paper are as follow: 
(1) present four new formulations for the FCMVRP, (2) 
compare the formulations both analytically and empirically, 
and (3) through extensive computational experiments, show 
that instances with maximum of 40 targets are within the 
computational reach of a branch-and-cut algorithm based on 
the best of the four formulations. 

The rest of the paper is organized as follows. Sec. Ill states 
the formal definition of the problem and introduces notations. 
In Sec. IV, we develop the four mixed integer linear pro¬ 
gramming formulations. The first two formulations are arc- 
based and the rest are node-based formulations i.e., decision 
variables for enforcing the fuel constraints are introduced for 
each edge and each target for the arc-based and the node- 
based formulations, respectively. The linear programming 
relaxations of the formulations are analytically compared in 
this section. In Sec. V, we present the computational results 
followed by conclusions and possible extensions. 

III. Problem definition 

Let T denote the set of targets {U,..., . Let D denote 

the set of depots or refueling stations {di,... ,dk}', each 
depot dk is equipped with a vehicle Vk. The FCMVRP is 
defined on a directed graph G = {V,E) where V = T U D 
and E is the set of edges joining any two vertices in V. 
We assume that G does not contain any self-loops. Each 
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edge {i,j) G E is associated with a non-negative cost 
required to travel from vertex i to vertex j and fij, the fuel 
spent by traveling from i to j. It is assumed that the cost of 
traveling from vertex i to vertex j is directly proportional to 
the fuel spent in traversing the edge {i,j) i.e., Cij = K ■ fij 
(cij and Cji may be different, but for the purpose of this 
paper, we assume Cij = cji). It is also assumed that travel 
costs satisfy the triangle inequality i.e., for every i,j,kG V, 
Cij -I- Cjk > Cik. Furthermore, let F denote the fuel capacity 
of all the vehicles. The FCMVRP consists of finding a route 
for each vehicle such that the vehicle Vk starts and ends its 
route at its depot dk, each target is visited at least once by 
some vehicle, the fuel required by any vehicle to travel any 
segment of the route which joins two consecutive depots in 
the route must be at most equal to F, and the sum of the 
cost of all the edges present in the routes is a minimum. 

IV. Mathematical formulations 

This section presents four formulations for the FCMVRP. 
The first two formulations are arc based, and the remaining 
formulations are node based. The arc based and edge based 
formulations have additional decision variables for each edge 
and vertex respectively, to impose the fuel constraints. For 
any given formulation F, let denote its linear program¬ 
ming relaxation obtained by allowing the integer variables 
to take continuous values within the lower and upper integer 
bounds, and opt(J^) denote the cost of its optimal solution. 

A. Arc-based formulations 

We first present an arc based formulation Fi for the 
FCMVRP, inspired by the models for standard routing prob¬ 
lems [24], [17]. Each edge {i,j) G E is associated with a 
variable Xij, which equals 1 if the edge {i,j) is traversed 
by the vehicle, and 0 otherwise. Also, associated with each 
edge {i,j) is a flow variable Zij which denotes the total fuel 
consumed by any vehicle as it starts from a depot to the 
vertex j, when the predecessor of j is i. Using the above 
variables, the formulation Fi is given as follows: 

(J^i) Minimize E C-ij 

{i.3)&E 


subject to: 



^ ^ ^di ~ ^ ^ ^id V d G D, 

i&V i&V 


(1) 

Xij = 1 and Xji = 1 

i&V i&V 

Vj G T, 

(2) 

E/ ~ E/ ~ E/ 

3&V j£V jev 

Vi G T, 

(3) 

0 < Zij < Fxij V {i,j) G E, 


(4) 

^di — fdidtdi Vi G T, d G D, 

and 

(5) 


XijG{Q,l) 'i{i,j)GE. 

In the above formulation the Eqs. (1) - (2) impose the degree 
constraints on the depots and the targets. The constraints in 
(3) are the connectivity constraints; they eliminate sub tours 
of the targets. (4) and (5) together impose 0 < Zij < F and 
they ensure that the fuel consumed by the vehicle to travel 
up to a depot does not exceed the fuel capacity F. Einally, 
the constraints in Eqs. (6) impose the binary restrictions on 
the variables. 

Now, we present another arc-based formulation F 2 which 
is a strengthened version of Fi. The following proposition 
is a modified version of the Prop. 1 presented in [17] for the 
distance constrained vehicle routing problem; it strengthens 
the bounds given by the constraints in (4). 

Proposition 1. The inequalities in (4) can be strengthened 


as follows: 




Vj G T, {i,j) G E 

(7) 

Zid fi Fxid yi G T and d G D 

(8) 


yi G T, {i,j) G E 

(9) 


where, U = mindg^, fid and Si = min^gc fdi. 

Proof. When j is a depot, the constraints in (8) and (4) 
coincide. We now discuss the case when both i and j are 
targets. When Xij = 1, any vehicle that traverses this edge 
{i,j) consumes at least (si + /ij) amount of fuel. As a result, 
the constraint in (9) strengthens the lower bound of Zij in 
(4). Similarly, the total fuel consumed by any vehicle that 
traverses the edge {i,j) cannot be greater that (F—tj), where 
tj is the minimum amount of fuel required by any vehicle 
to reach a depot from target j. Therefore, the constraint in 
(7) strengthens the upper bound of Zij in (4). □ 

Hence, the second arc-based formulation is as follows: 
{F 2 ) Minimize E Cij Xij 

subject to: (1) - (3), (5) - (6), and (7) - (9). 

Corollary 1. opt(J'^) > opt(J^j^). □ 

B. Node-based formulations 

In this section, we present a node-based formulation for 
the ECMVRP based on the models for the distance con¬ 
strained VRP in [25], [16]. Eor the node based formulation, 
apart from the binary variable Xij for each edge {i,j) G E, 
we have an auxiliary variable Ui for each vertex i, that 
indicates the amount of fuel spent by a vehicle when it 
reaches the vertex i. We assume = 0 as the vehicles 
are refueled to their capacity when they reach a depot. In 
addition, we will also use the following two parameters: 
ti = mindgn fid and s* = min^gc fdi for every vertex 



i ^ V. For any d G D, td = 0 and = 0. Using the 
above notations, the formulation T 3 is given as follows: 

( 7 ^ 3 ) Minimize ^ CijXij 


(i.peE 



subject to: (1), (2), and (6), 



Ui — Uj + MXij < M — fij 

yiGV,jGT, 

(10) 

tti ^ V ^ ^ ifdi 
dGD 

Vi e T, and 

(11) 

Ui<F-U - ^{fid - U] 

dGD 

\Xid Vt e T. 

(12) 


The constraint in Eq. (10) serves both as sub tour elimination 
and fuel constraints. It eliminates sub tours of the targets and 
ensures any route that starts and ends at a depot consumes 
at most F amount of fuel. This can be easily observed by 
aggregating the constraints for any sub tour of the targets and 
for any route starting and ending at a depot [25]. The value 
of M in the constraint is given by M = max(j ,j)eE{F - 
Sj — ti + fij}. The constraints in Eqs. (11) and (12) specify 
the upper and lower bounds on Ui, for every vertex i. The 
following proposition strengthens the fuel constraints and the 
bounds on Ui. 

Proposition 2. The inequalities in (10), (11), and (12) can 
be strengthened as follows: 

Ui - Uj + Mxij + {M - fij - fji)xji < M - fij 


y^,jyT, 

(13) 

Ui ^ ^ (.Sj f ji) X ji yi G T, 

(14) 

jeu 


Ui < F — ''^^{tj + fij)xij Vi G T, and 

(15) 

j&v 


Ui<F-U-^{F-ti- fdi)xdi Vi e T. 

(16) 


deD 

where, xu = 0 and Xij = 0 whenever Si + fij + tj > F. 


Proof The constraint in Eq. (13) can be obtained by lifting 
the variable Xji in Eq. (10). A constraint is said to be “valid” 
if it does not remove any feasible solution to the ECMVRP. 
We compute the value of the coefficient a that makes the 
following constraint valid: 

Ui — Uj + Mxij + axji < M — fij. 

The equation is valid when Xji = 0, as it reduces to (10). 
When Xji = 1, we have Xij = 0 and Uj + fji = Ui. Hence, 
the best value of a that makes the equation valid is given by 

M - 

Similarly, Eq. (14) can be obtained by lifting every Xji 
variable for j gT m any order. We will illustrate the lifting 


procedure for one of the Xji variables. This involves com¬ 
puting the coefficient a that makes the following constraint 
valid: 

^ -Si “f ^ ^ {fdi ^i^tXdi CXXji. 
dGD 

The above equation is valid when Xji = 0, and when Xji = 1, 
we have Xdi = 0 and a < Ui — Si. The best value of a that 
does not remove any feasible ECMVRP is hence given by 
Sj + fji — Si. Similarly, the coefficients of the other Xji 
variables can be computed. The resulting constraint is given 
by 

tti ^ -f ^ F fji ’^fjtCji Vi G U. 
j&v 

In the above equation, Sj = 0 for j G D. The above equation 
reduces to Eq. (14) due to the degree constraints in (2). The 
constraints in Eq. (15) are similarly obtained from (12) by 
lifting the Xij variable for every j e T. The proof is omitted 
as it is similar to the previous ones in the proposition. The 
constraints in Eq. (16) are valid bounding constraints for the 
ECMVRP when the target i is the first target that is visited 
by any vehicle as it leaves the depot. In this case, the Eq. 
(12) reduces to < F — ti. We further strengthen this 
constraint by lifting the variable Xdi for every d G D. The 
lifting coefficient a for Xdi takes the value —{F — ti — fdf) 
and the resulting constraint is given by Eq. (16). □ 

Hence, the second node-based formulation is as follows: 

{Ff) Minimize E ^ij ^ij 

subject to: (1), (2), (6), and (13) - (16). 

Corollary 2. opt(J'.f') > opt(J^I'). □ 

V. Computational results 
In this section, we discuss the computational performance 
of the four formulations presented in the previous section. 
The mixed integer linear programs were implemented in 
Java, using the traditional branch-and-cut framework of 
CPLEX version 12.4. All the simulations were performed 
on a Dell Precision T5500 workstation (Intel Xeon E5630 
processor @2.53 GHz, 12 GB RAM). The computation 
times reported are expressed in seconds, and we imposed a 
time limit of 3,600 seconds for each run of the algorithm. 
The performance of the algorithm was tested with randomly 
generated test instances. 

Instance generation 

The problem instances were randomly generated in a 
square grid of size [100,100] with 5 fixed depot locations. 
The number of targets varies from 10 to 40 in steps of 
five, while their locations were uniformly distributed in 
the square grid; for each \T\ G {10,15,20,25,30,25,40}, 
we generate five random instances. Each depot contains a 
vehicle. The travel costs and the fuel consumed to travel 



between any pair of vertices are assumed to be equal to 
the euclidean distances between the pair. For each of these 
problems, we generate four possible fuel capacities F as 
a function of the the distance to the farthest target from 
any depot A. The fuel capacity F of the vehicles gets the 
values 2.25A, 2.5A, 2.75A and 3A. In total, we generate 140 
instances and run the branch-and-cut algorithm for all the 
formulations. 

Tables I and II, and Fig. 2-3 summarize the computational 
behavior of the algorithms for all the 140 instances. The 
following nomenclature is used throughout the rest of the 
paper: 

instance number; 

opt(J'/^): linear programming relaxation solution for formu¬ 
lation 

n: instance size i.e., number of targets in the instance; 
%-LB: percentage LB/opt, where LB is the objective value 
of the linear programming relaxation computed at the root 
node of the branch and bound tree and opt is the cost of the 
optimal solution to the instance; 
total: total number of test instances of a given size; 
succ: number of instances for which optimal solutions were 
computed within a time limit of 3,600 seconds. 

Table I compares the cost of the linear programming (LP) 
relaxations of the four formulations presented in Sec. IV for 
the 40 target instances. The results in table I provide an 
empirical comparison of the formulations presented in IV; 
the observed behavior is expected because the formulations 
F 2 and F 4 are strengthened versions of Fi and F 3 , respec¬ 
tively (see corollaries 1 and 2). As for the LP relaxations of 
formulations F 2 and F 4 , it is difficult to conclude that one 
is better than the other since F 4 produces better relaxation 
values than F 2 only for 60% of the instances. Hence, the 
rest of the computational results compares the formulations 
F2 and F4. 

Table II shows the number of instances of different sizes 
solved to optimality by the formulations F 2 and F 4 within 
the time limit of 3600 seconds. The plot in Fig. 2 shows the 
average time taken by the two formulations to compute the 
optimal solution to the FCMVRP. The table II and Fig. 2 
indicate that the arc-based formulation F 2 outperforms the 
node-based formulation F 4 for the larger instances. For the 
smaller sized instances, it is difficult to differentiate between 
the two formulations. The plot in Fig. 3 shows the percentage 
LB/opt for both the formulations (LB is the objective value 
of the linear programming relaxation computed at the root 
node of the branch and bound tree and opt is the cost of 
the optimal solution to the instance; for the instances not 
solved to optimality, opt represents the cost of the best 
feasible solution obtained at the end of 3,600 seconds). We 
observe that the %LB is consistently better for formulation 
J^ 2 - This plot also provides empirical evidence to the claim 


TABLE I 

Cost of the LP relaxation for the 40 target instances. 


# 

optjj'f) 

opt(J'l') 

OPV-^f) 


1 

496.42 

509.24 

426.17 

518.00 

2 

487.31 

496.39 

426.17 

518.00 

3 

480.55 

487.40 

426.17 

518.00 

4 

475.23 

480.33 

426.17 

518.00 

5 

444.35 

458.01 

389.08 

434.00 

6 

435.45 

445.70 

389.08 

434.00 

7 

428.44 

436.47 

389.08 

434.00 

8 

423.06 

429.97 

389.08 

434.00 

9 

396.10 

403.96 

367.11 

452.00 

10 

392.87 

398.72 

367.11 

452.00 

11 

390.42 

394.66 

367.11 

452.00 

12 

388.40 

391.85 

367.11 

452.00 

13 

481.22 

493.64 

427.04 

461.00 

14 

469.76 

479.81 

427.04 

461.00 

15 

461.16 

469.20 

427.04 

461.00 

16 

454.80 

461.47 

427.04 

461.00 

17 

503.19 

516.58 

461.07 

523.00 

18 

494.98 

504.84 

461.07 

523.00 

19 

489.64 

496.31 

461.07 

523.00 

20 

485.92 

489.99 

461.07 

523.00 


TABLE II 

Comparison of formulations J2 and J4. 


n 


F2 

Fi 

total 

succ 

succ 

10 

20 

20 

20 

15 

20 

20 

20 

20 

20 

20 

20 

25 

20 

20 

14 

30 

20 

20 

5 

35 

20 

20 

15 

40 

20 

19 

1 


that the arc based formulation F 2 outperforms the node based 
formulation F 4 . 

VI. Conclusions and future work 

In this paper, we have presented four different mixed 
integer linear programming formulations for the multiple 
depot fuel constrained multiple vehicle routing problem. The 
problem arises frequently in the context of path planning for 
UAVs, green vehicle routing and routing electric vehicles. 
The formulations have been compared both analytically and 
empirically, and it is observed that a strengthened arc based 
formulation (F 2 ) performs the best in terms of computing 
optimal solutions to the problem. Computational experiments 
on a large number of test instances corroborates this observa¬ 
tion. Future work can be directed towards developing similar 
mixed integer linear programming formulations and branch- 
and-cut algorithms to solve a heterogeneous variant of the 
problem i.e., with vehicles having different fuel capacities. 
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